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STABILITY ESTIMATES FOR THE REGULARIZED INVERSION OF 
THE TRUNCATED HILBERT TRANSFORM 

RIM A ALAIFARI 1 , MICHEL DEFRISE 2 , AND ALEXANDER KATSEVICH 3 


Abstract. In limited data computerized tomography, the 2D or 3D problem can be 
reduced to a family of ID problems using the differentiated backprojection (DBP) method. 
Each ID problem consists of recovering a compactly supported function / £ L 2 (J r ), where 
T is a finite interval, from its partial Hilbert transform data. When the Hilbert transform 
is measured on a finite interval Q that only overlaps but does not cover T this inversion 
problem is known to be severely ill-posed [I]. 

In this paper, we study the reconstruction of / restricted to the overlap region T n Q■ 
We show that with this restriction and by assuming prior knowledge on the L 2 norm or 
on the variation of /, better stability with Holder continuity (typical for mildly ill-posed 
problems) can be obtained. 


1. Introduction 

This paper presents new results on the stability of the inversion of the Hilbert transform 
with limited data. The main application concerns tomographic reconstruction from trun¬ 
cated projections, where some 2D and 3D problems can be reduced to the inversion of the 
Hilbert transform on a family of line segments defined by the geometry of the scanners. The 
reduction of the 2D or 3D tomography problem to a ID Hilbert transform inversion prob¬ 
lem is achieved by backprojecting a derivative of the projection data over an angular range 
of 180 degrees. This operation based on a result by Gelfand and Graev [Sj is referred to 
as the differentiated backprojection (DBP). In the 2000’s, the DBP triggered a remarkable 
evolution in tomography by allowing accurate reconstruction of regions of interest (ROI) 
from truncated projection data sets mmm, which were previously believed to exclude 
any accurate reconstruction. A second method for accurate ROI reconstruction from lim¬ 
ited data, the virtual fan-beam algorithm (VFB), was introduced in 2002; it is based on 
a different mathematical property of the Radon transform and will not be studied in this 
paper. The reader is referred to [6| for a review on ROI reconstruction in tomography and 
the link between that problem and the Hilbert transform. 

We work with a line integral model of the tomographic data, assuming continuous (infin¬ 
itely fine) sampling of the 2D or 3D x-ray transform of the unknown function f(x) within 
some subset of the space of lines in M 2 or M 3 . In this context ’’accurate” reconstruction 
means that the data uniquely and stably determine f(x) within a ROI, under appropriate 
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assumptions. These assumptions specify the functional spaces to which the object and the 
data belong, the nature and magnitude of the measurement noise, and typically also involve 
prior constraints on the object. Loosely speaking we call a reconstruction ’’stable” with 
Holder continuity if the norm of the error on the estimated f{x) (possibly within some 
limited ROI) can be bounded by some positive power of the norm of the measurement 
error [5j. In contrast one speaks of logarithmic continuity when the error bound decreases 
only with the logarithm of the noise. The noise propagation properties of the DBP are 
well understood, since they are similar to those of the usual filtered-backprojection (FBP) 
algorithm. This can be seen by observing that the two procedures involve the same back- 
projection, and that the ramp filter \v\ (for FBP) and the derivative filter 2ixiv (for DBP) 
similarly amplify high frequencies. Thus, the major question left is to analyze the stability 
of the inversion of the Hilbert transform with limited data. The rest of this paper focuses 
on that one-dimensional inverse problem. 

To date, much more is known on conditions guaranteeing uniqueness of the solution than 
on the stability in the presence of noise. Stability is easy to analyze when a closed form 
analytic inversion formula is known, such as the FBP algorithm for the Radon transform 
with complete data. With the DBP however, a closed form inversion exists only when 
the limited data allow calculating the Hilbert transform of f(x) on a line segment Q that 
contains the support T of fix) along that line. In this case, inversion is based on the finite 
inverse Hilbert transform and has very good stability (all singular values but one are equal 
to 1 in a weighted L 2 space [IT])- Unfortunately no closed form inversion is known for the 
two other configurations: the interior problem (Q C T) and the truncated problem where 
Q partially overlaps T. 

It has long been known m that unique inversion is not possible for the interior prob¬ 
lem, even within the ROI Q where the data are known. However, uniqueness holds when 
additional information is available. One such instance is the case where f(x) is known in 
some segment K, C Q mmm- Other results show that uniqueness is restored when f(x) 
is assumed to belong to some function space (e.g. f(x) is piecewise constant, a polynomial 
function or a generalized spline) |18J [19]. Numerical tests with discretized models of these 
problems suggest good stability but these results may depend on the way the problem is 
discretized. To our knowledge explicit stability estimates have only been obtained when 
f(x) is a polynomial function in the interior region [11 ]. 

For the truncated problem uniqueness of the solution follows from simple analyticity 
arguments [8] . However, the stability bound obtained in [8] is difficult to interpret since it 
relates the reconstruction error to an upper bound on the error on an intermediate function 
and not to the error on the Hilbert transform. Also, the results in |8j do not pertain to 
a specific regularization method. Recently, a more explicit stability estimate has been 
obtained in [3(j. It guarantees logarithmic continuity for the reconstruction on the entire 
support T of fix') if prior knowledge on the total variation of / is assumed. 

This paper exploits recent results on the asymptotics of the singular value decomposition 
(SVD) of the truncated Hilbert transform with overlap [lj. With these, explicit l? stability 
estimates are obtained. They allow to guarantee stable inversion of the corresponding 
inverse problem when an a priori bound on the L 2 norm of the solution is assumed. We 
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give specific examples for the truncated SVD and for Tikhonov’s method. By seeking to 
reconstruct only within a region of interest , Holder continuity is obtained. A similar result 
is shown for prior knowledge on the total variation of f(x). 

Remark. Stability estimates allow assessing the ill-posedness of an inverse problem 
independently of the algorithm and independently of the discrete sampling of the data 
and solution. The inversion of the truncated Hilbert transform is a very specific problem, 
which is severely or mildly ill-posed depending on which region of interest is reconstructed. 
The stability estimates in this paper are derived in a deterministic framework, they give 
worst case error bounds, which are pessimistic and for that reason do not provide reliable 
absolute values for the regularization parameter. However, our results for the Hilbert 
transform determine the Holder power of the reconstruction accuracy in terms of known 
constants depending on the geometry of the problem (i.e. on the intervals J- and Q). This 
may provide guidance to adapt the regularization parameter to varying geometry or data 
error, once it has been empirically estimated for one case. Future work will investigate 
application to image reconstruction from incomplete CT data. 

2. The truncated Hilbert transform 

Consider the truncated Hilbert transform H? : L 2 (J r ) -A L 2 (Q), with T = (02,04), 
Q = (01,03), and ai < 02 < 03 < 04: 

(2.1) ( H T f)(x ) = -p.v. f ^ y ^ dy , x £ (03,03) 

Ja 2 V - x 

We define the normalized singular function pair u n £ L 2 (J 7 ), v n £ L 2 (Q) such that Htu h = 
a n v n for n £ Z. The singular values a n are ordered as usual by decreasing values, and 
the spectrum has two accumulation points, lim n ^._ 0O a n = 1 and lim^^oo a n = 0 . We will 
make use of the following properties of the singular system of Ht- 

- For large n > 0 the asymptotic behavior of the singular values is given by 

(2.2) a n = 2 e-™ K+/K - (1 + 0(n” 1/2+c )), n -A 00, 


for some small £ > 0 and with 



Here, P{x) = (x — a\)[x — a2){x — 0,3)(x — 04). 

- For large negative n the asymptotic behavior of the singular values is 

( 2 . 4 ) a n = (1 - 2 e~ 2 ^ nK - /K +)(l + 0 (|n|- 1/2+? )), n -A -00. 

The singular function v n is bounded at ai and 03 and has a logarithmic singularity 
at a/2- For large n > 0 it is an oscillatory function in the interval (04,02) and a 
monotonically decreasing function on (02,03). 

The function u n is bounded at 02 and 04 and has a logarithmic singularity at 03. 
For large n > 0 it is a monotonically increasing function in the interval (02, 03) and 
an oscillatory function on (03,04). 
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For sufficiently small /x > 0, the norm of u n over the segment (02, 03 — /x) has the 
following asymptotic behavior 

1/2 


( 2 . 5 ) 


( 2 . 6 ) 


r a 3-m 


dx |u n (x)|' 


'Q2 


1 


where 




r a 3 


dt 


/mr 


2 vr 


^-^"(l + Ofa- 1 / 2 -*)), 


v/M 


■ (1 + 0(h)). 


“3—/x VTO T 3 ^) 

Except for the last property which is proven in the appendix (see Lemma [I]), the proofs of 
all the above properties can be found in mm- However, the nronotonicity of the singular 
functions u n is shown more explicitly in the appendix (see Lemma [2]). 

3. Inversion of the truncated Hilbert transform: regularization with an 

L 2 penalty. 

We consider the following inverse problem. Let f ex E L 2 (T") be an object of interest 
and let g ex = B^fex £ L 2 (Q) be the corresponding noise-free data, with Ht the truncated 
Hilbert transform ( 2 . 1 ). Given a noisy measurement g such that \\g — gex\\L 2 (g) — & for 
some noise level 5 > 0 , we seek an estimate of f ex on some interval (<22, <23 — /x), with some 
small g > 0. 

The rationale to restrict the reconstruction to a limited interval is that there is no hope 
of obtaining a Holder stability estimate outside the segment Q = (ai, 03) where the Hilbert 
transform is known. Following a similar intuitive analysis by Natterer [ 2 ] for the exterior 
problem of tomography, this can be seen by considering an object f ex which is supported 
in (03,124) and has a discontinuity somewhere in that interval. For such an object, one 
sees from (2.1) that Hxfex is C°° and therefore an inverse operator should map a C°° 


function onto a discontinuous function, which corresponds to a severely ill-posed problem. 
This severe ill-posedness is also expected from the exponential decay of the singular values, 
equation ( 2 . 2 ). Thus, unless very restrictive prior knowledge is assumed the quest for a 


stable reconstruction must be restricted to an interval such as (<22,03-/2). By varying the 
parameter /x one can study the degradation of the stability as one gets closer to the limit 
of the stably recoverable region (02,03). We denote by %/x the characteristic function of 
that ”region-of-interest” interval (02,03 — /x). 

We assume the following properties of the SVD: 

- PI: 0 < a n < 1 . 

- P2: There are some Nq e N and positive real A, a such that cr n > Ae~ an for 
n > N 0 . 

- P 3 : There are some N fl e N and positive real such that \\Xn u n\\L' 2 (T) — 

for n > N\ 1 . 


The constants A and a can be obtained from ( 2 . 2 ): 


(3.1) 


A < 2 , 


a = ttI \ + / , 
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where A = 2 corresponds to the leading term of the asymptotic behavior of the singular 
values. Hence, any value of A smaller than 2 provides an upper bound for sufficiently large 
n. We note that the index Nf and the constants B tl , ( 3 ^ all depend on the choice of g. In 

_ i 

particular, for P 3 to hold, = 0 (g 1+2C ); so the smaller g is chosen to be, the larger N^ 
has to be taken. Throughout this paper, we assume g > 0 to be arbitrarily small but fixed. 
The results obtained in the following sections are no longer valid as g —> 0 . This is due 
to the fact that Wx^oUnW^^ only decays as 0 (l/n) (see equation (8.6) in (lj). In what 
follows, we select N tl such that Nf > Nq, so that P 2 holds for all n > N^. Property P 3 
follows from equation ( 2 . 5 ), this asymptotic relation can be changed into an upper bound 
by taking for instance 

1 


vT 


tt 


3 . 1 . An algorithm independent stability bound. We assume as regularizing prior 
knowledge on the solution an upper bound E > 0 on the norm of the object, \\fex\\L' 2 (T) — 
E. The stability bound that we present follows the approach of Miller m • This formulation 
is independent of a specific algorithm, but instead yields an upper bound on the L 2 distance 
(within the interval (02,03 — g)) between any pair of solutions that are compatible both 
with the data and with the prior knowledge. 

We define the set of admissible solutions as 


5 


{/ G L\E) 


II H T f ~ g\\ L z(G) < 5 and \\f\\ L 2( T) < E } . 


Proposition 1 . Let Erf ex = 9ex and assume \\fex\\L 2 (T) < E. We consider the recon¬ 
struction of f ex from a noisy measurement g for which \\g — g e x\\L 2 (g) < d is given. Then, 
any algorithm that for given 5 , E and g computes a solution in S, is a regidarization method 
for the reconstruction on (02,03 — g). More precisely, given any two admissible solutions 
f a , fb G S, one has for sufficiently small 6 the bound 


( 3 . 2 ) 


2 de aN » 

\\x,*(fa-f b )\\L*(T)<—j- + 2 EB li 


( a 

AV^E J \ (a - fd/x) yj e 2 ^ - 1 


where V fl is a constant, which only depends on a and ( 3 ^. 


Proof. Consider two admissible solutions f a , fb G <S. The corresponding Hilbert transforms 
9a = H T f a and g b = H T f b satisfy ||g a - 9 b\\L 2 {g) < ha - g\\L 2 {g) + II 9 b ~ g\\L 2 (g) < 25 . 
Similarly ||/ Q — f b \\L 2 (r) — 2 E. We can rewrite f a — f b using the SVD of Ht and splitting 
the obtained expansion into three terms as follows: 

Alp 1 AT j v 

(3.3) fa~ fb = ^ ' ( 9a 9b,’an) ^nT ^ ' ( 9a 9b, Vn) T ^ ' ( fa~fb,'a n )u n . 

n=—co n=N, + l n^+1 


Here, the SVD series is split at and at some - so far arbitrary - index N and we used 
{9a, Vn) = { H T fa,v n ) = cr n {f a ,u n ) and idem for f b . Equation (fT 3 |) holds for any N > N^. 
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We seek to derive an upper bound on the error in the interval of interest in the form of 
( 3 - 4 ) II X/iifa — fb)\\L 2 (F) < h + h + h 


where I\, I2,13 correspond to the norms of the three terms in ( 3 . 3 ) restricted by \n- 
The error due to the large singular values is bounded as 


l N » 1 

if = IIXm E ( 9a ~ 9b > Vn ) — Un Wh(T) ^ II E ( 9a ~ 9h > Vn ^> ~ 

(J n (J n 


Ur 


Il 2 (^) 


1 ^ 

( 3 . 5 ) <~ 2 ~ E \( 9 a- 9 biV n )\ 2 < 

^Nn rt.=—n 


" a k 


n=—o o J 

The contribution of the intermediate terms corresponding to N^ < n < N to the error 
is bounded as 

2 


r 2 _ 
1 2 ~ 


N 1 j N 1 

E {9a-gb,V n ) X^Un\\ 2 L 2 {F) < { Y \{9a-9b,Vn)\ 11 XfJ.U n \ \ L 2 ( F ) 

O'n. I .. U .> 


n=Nfj +1 
TV 


n=Nu+l 


TV A 

< E I (9a - 9b, V n >) | 2 ^2 -^\\X^Un\\ 2 L 2 ^) < (26) 2 Y T 11X/R^ 11Z, 2 (.F) 


< 


(2 5 ) 2 Bl 

A 2 


TV 

E 

n=TV„ +1 


o 2(a-^)n _ ( 2J ) 2 ^ e 2(a P^ N -e 2 ( a 


A 2 


1 — g 2(a ftn) 


( 3 . 6 ) 


< 


( 2 < 5) 2 B 2 e^a-MN 


A 2 1 — e -2(a-^) ' 


where we successively used the triangle inequality, Schwarz inequality, and majorized by 
extending the sum over n' to all Z. Recall from the definitions of these constants that 
a > / 5 fj > 0 . The last inequality is not strictly needed but simplifies subsequent derivations, 
and it is expected to have a limited impact since for small 5 , we will find that the optimal 
cut-off satisfies N 3 > N^. 

Finally, 


^3 — II /* 1 (fa fbi u n) < S 'Y. l(/a fbi u n)\ HX/i^n lli 2 ^) Z 

n=N-\-l vn=7V+l J 

00 00 00 

— El \(fa~fbi u ri) | IIX^nH^pr) < \\fa — fb\\L 2 (T) El IIX/i w n |lx,2(jr) 

7i'=./V+l n=iV+l n=iV+l 

00 


( 3 . 7 ) < ( 2 E) 2 Bl Y e“ 2 ^ n = { 2 E B^e~^ N /ke 2 ^ - l}" 


n=N +1 
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Inserting equations (3.5), (3.6), (3.7) into (3.3) leads to the following upper bound on 
the L 2 distance between f a and fb on the region-of-interest interval ( a 2,03 — /r): 


(3.8) 


\\Xn(fa — fb)\\L 2 (T) < 


2d 25B U 


e (a-MN 


a N u 


A 


+ 


2 EB^e~^ N 


\J 1 — e 2(o: /3 m ) yj \ 


A quasi-optimal splitting index N (d) is obtained by treating N as a real variable^] and by 
minimizing the convex function in the RHS of (3.8) w.r.t. N: 


(3.9) 

with the constant 


N(d) = - log 
a 


EAV„ 


V = _^_ 

11 (a - /3„) 


(1 


_ _-2(a-0, 




1/2 


(e 2 ^ — l) 


Inserting IV(d) into (3.8) finally yields (3.2). This error bound is valid when N(5) > Nn, 


a condition which is always satisfied for sufficiently small noise level d. 

To conclude the proof, consider any algorithm that produces for each d, E and g a 
solution f a E S. The bound (3.2) can be applied with fb = f e x because f ex E S, and since 
the right-hand side in ( |3.2[ ) tends to zero as d —> 0, this algorithm is a regularizing method 
for the reconstruction on ( 02 , 03 — //). □ 

Remarks. 


For small d the second term in (3.2) eventually dominates. An intuitive interpreta¬ 


tion of the stability estimate (3.2) is that the number of significant digits that can 


be reliably recovered when estimating the solution is roughly equal to a fraction 
/3n/a of the number of significant digits in the measured data. 


The first sum in (3.3) starts at —00 because the spectrum of Ht has two accu¬ 
mulation points. We remark that the singular components up to n —> —00 can be 
stably recovered from the data because a n —> 1 as n —> — 00 . This holds only for 


the continuous-continuous model of the inverse problem, which is analyzed in this 
paper. In practice, only sampled data are available. Since the functions u n are 
increasingly oscillating functions on the interval (< 12 , 0 - 3 ) when —> —00 lD( 2 j, one 
would in practice start the SVD expansion at a large negative N m i n < 0. This lower 
SVD cut-off would have an effect similar to a cut-off on the high spatial frequencies. 
Setting g = fb = 0 leads to another equivalent formulation of stability: if / E L 2 {E) 
is such that \\Ht fW^^g) < 5 and ||/|| L 2 ^ < E, then 


(3.10) 


Se aNfi 

IIXm/IIt 2 (-F) < —^—h eb ^ 


AV^E 


h lo¬ 


ot 


(a - j3n) Ve 2 ^ - 1 


^We hereby neglect the fact that N only takes integer values. This hardly changes the error estimate 
when S is small. 
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The following corollaries of Proposition [I] give explicit error bounds for two specific 
regularization algorithms. The proofs follow the same line as the proof of Proposition [l 
and are therefore omitted. 

Corollary 1 (Truncated SVD). The truncated SVD estimate 

1 


(3.11) 


N{5) 

fN(5) = ^2 


Vn) V"n 


with the cut-off index N(5) defined by (3.9) guarantees regularization in the sense that 
(3.12) 

Specifically, if N ( 6 ) > N 


Jim II Xti{fN(5) ~ fex)\\ L 2(T) = 0. 


(3.13) WxMn - fex)\\ L Z(T) < 


5e aN » 


+ EB Ll 


AV^E 


Pn/ a 


a 


(a - fiA Ve 2 /V - 1 ) 


Note that the quasi-optimal number (3.9) of terms in the truncated SVD (3.11) increases 
only very slowly with decreasing noise level 6 , which is a typical behavior of severely ill- 
posed inverse problems. The corresponding lower bound on the singular values that are 
included in the truncated SVD solution estimate, 


(3.14) 


is proportional to the noise to signal ratio 5/E. This behavior of the cut-off singular value 
is also obtained for inverse problems with compact operators using Miller’s method 0113!, 
0 p. 258], 

Corollary 2 (Tikhonov regularization). The Tikhonov estimate 


( 3 . 15 ) f v = arg min <Z> v (f) with $„(/) = ||iZr/ - 9\\b{g) + v\\f\\b(,r) 

with r] = r](5) > 0 such that rj(5) —> 0 and that 5 2 /r](5) remains bounded as 5 
guarantees regularization in the sense that 

J™ \\Xn(fv(S) ~ fex)\\L*(F) = 0- 

In particular the choice rj = 5 2 / E 2 leads to 

II Xuifv ~ fex)\\L 2 (F) < (1 + V2 )——— 

( 3 . 16 ) +(1 + V2 )EB„ 


0 , 


AV^E 


/W Q 


a 


(a - y/e 2 P» - 1 


provided 5 is small enough such that N(5) > IV^ in (3.9). 
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The exponential decay of the singular values of the truncated Hilbert transform (see 


equation (2.2)) is typical of severely ill-posed inverse problems. Equations (3.2 3.13, 3.16), 
however, give a Holder continuity expected for mildly ill-posed problems. This mild ill- 
posedness is due to the fact that the error is only calculated over a restricted interval 
( 02,03 — n), which excludes a neighborhood of the limit 03 of the segment where the 
Hilbert transform is known. The Holder error bound does not hold in the limit /r —> 0; 
the squared norm \\ Xu = o u n \\ 2 L 2 ^ J r ' ) decays as 0 ( l / n ) (equation ( 8 . 6 ) in [ 1 ]), too slowly to 
balance the exponential decay of a n . 

The error bounds obtained have the same dependence on /“. This behavior is de¬ 
termined essentially by the exponential decay rates a and (3^, but also by the assumed 
regularizing prior knowledge ||/||l 2 (.f) < F. Figure [l] illustrates the change in the rate 
j/W" as the amount of overlap between the two intervals T and Q is varied. In the follow¬ 
ing section, we seek to find the Holder exponent using a non-quadratic regularizing penalty, 
or more precisely, a total variation (TV) penalty. It turns out that this yields the same 
dependence on <5^/ a . 



Figure 1. The Holder power /3^/a is plotted for 0 < 03 < 1 , and for fixed aj = 
— 1 , 0,2 — 0 , 0,4 — 1 to illustrate the dependence of the convergence rate of the regularized 
solution on the size of the overlap region ( 02 , 03) relative to the function support ( 02 , (24) = 
(0,1). The green, red and blue curves correspond to u = 0.25 • < 23 , 0.1 ■ <23 and 0.01 ■ (23. 


4. Inversion of the truncated Hilbert transform: Stability through a TV 

PENALTY 

In this section, we drop the condition ll/exlli, 2 ^) < E and instead consider regularization 
of the inversion problem by assuming prior knowledge on the variation of f ex in the form 
of the upper bound \fex\rv < ft for some n > 0. Additionally, we assume that f ex vanishes 
at 02 and 04. This latter assumption is not too restrictive; for if f ex does not vanish at 
the boundaries, we can always artificially enlarge the interval T = (02,04) such that the 
support of f ex is a strict subset of it. 

We define the set of admissible solutions as 

(4.1) S = {/ G BV(F) : || H T f - fj\\ L 2 (g) < <5, \f\rv < « and / vanishes at a 2 ,a 4 }. 
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Proposition 2. Let Hxfex = hex and assume \f ex \TV < K and f e x{a,2) = fex(a 4 ) = 0. We 
consider the reconstruction of f ex from a noisy measurement g for which \\g — g ex \\L 2 (g) < £ 
is given. Then, any algorithm that for given 5, k and g computes a solution in S, is a 
regularization method for the reconstruction on ( 02 , a$ — p). More precisely, given any two 
admissible solutions f a , fb £ <S and for sufficiently small 5, one has the following bound: 


( 4 - 2 ) II Xuifa ~ /&)||l 2 (.f) < 2 <L 4 1 e aNy ‘ + — B^k « 


5 \~ 

u \ ot 


a 


AW^J (a-p)(eP»- 1 ) : 


Here, c is a constant only dependent on a\, ... ,04 and W M is a constant only dependent on 
a, j3n, c and N^. 


Proof. Let g a = Hj-f a and gi, = H^fb- Then, we can write the same identity as in (3.3) 
which is exact for any choice of N > 7V /t . Lemma 6.1 in [JJ implies the following decay rate 
for sufficiently large index n: there exists a constant c only dependent on a\,... , such 
that for n > Nq, 

1 if r \\ pp' \fa fb\TV 
if a - fb,Un}\ < C -. 

n 

Similarly to the derivations in the previous section, we can obtain an upper bound on 
||Xuifa ~ fb)\\ l 2 (t) < h + h + L 3 , where I\ and I 2 are bounded as in (3.5) and (3.6), 
respectively. An upper bound on the last term is obtained as follows: 

OO OO 

-^3 = || /* ' (fa ~ fbi u n)Xu u 'n\\L 2 (J r ) — ^ ^ \(fa ~ fbi u n)\\\Xn u n\\L 2 (J 7 ) 


n=N+l 


n=N +1 


<c\fa-fb\rv ^2 ~\\ Xu u n \\ L 2 ( T ) ^ c \fa - fbhvB^ ^ 

n = N + 1 n=/V+l 


n 


< 


2 cnB^e~^ N 
NJeP* - 1 ) ' 


The last inequality in the above is not sharp but simplifies subsequent derivations. With 
this we obtain 


\\Xfi(fa - /&)||l 2 (.f) < 


25 2 5B U 

+ 


e {a-P JN 


a N u 


+ 


2 cKBpe-h * 1 


A (1 - e -2(a-^))i/2 AT M (e^-l)‘ 


The value N that minimizes this upper bound is found to be 


_ ff, (1 - e~ 2(Q ~ /3 ' j) ) 1/2 c 
^ ct-fin N^ePn - 1 ) 


where 


This choice yields (4.2). This bound is valid if N(5) > AT, which can be reformulated as 
an upper bound on the ratio 5/k: 


(4.3) 


- < AW u e 


-aN. 
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To conclude the proof, consider any algorithm that for each 5, n and g produces a solution 
f a € S. The bou nd (|4.2[ ) can be applied with ft, = f ex because f ex £ S, and since the 
right-hand side in (4.2) tends to zero as <5 —> 0, this algorithm is a regularizing method for 
the reconstruction on ( 02 , <23 — g). □ 


Remark. The bounds in (4.2) together with (4.3) imply that the error within the ROI 
tends to zero as n —> 0 even for constant noise level 5. This is only due to the prior 
assumption that f ex vanishes at 02 and < 24 . Together with \f e x\rv = 0 this assumption 
implies that / = 0. Any stability estimate that only requires an upper bound on \fex\rv 
without assuming that f ex vanishes at the boundaries, will not tend to zero as a -A 0 . 


4.1. Reconstruction on the entire interval. So far, we have found stability estimates 
for the reconstruction on the overlap region (or region of interest) while stable reconstruc¬ 
tion outside of the region of interest has not been guaranteed. One might ask whether 
it is possible to guarantee stable reconstruction on all of T = ( 02 , 04 ). If the assumed 
prior knowledge is of the form ||/ex||L 2 (j') £ A, such an estimate cannot hold in general. 
This can be seen for example by taking the sequence of singular functions u n for which 
IM|l 2 (.f) = 1 while \\H T Un\\Li{Q) 0. 

In [3] it has been shown that stability on all of T = ( 02 , 04 ) can be guaranteed if the 
variation of the solution is assumed to be bounded, i.e. \f\rv < «, for some n > 0. The 
estimate obtained is of logarithmic continuity (as opposed to Holder continuity). This is 
typical for severely ill-posed problems. 

We seek to derive such an estimate employing the methods applied in the previous 
sections. More precisely, we find an upper bound on \\f a — fb\\ i?(t) for fa,fb£<S where S 
is the set of admissible solutions in (14. lh. For this, we write 


II fa ~ fb\\L 2 (T) < h + h, 


where 


N 


— II 'y ^ (da 9b: v m) UuWl 2 ^) 
(J n 


and N > Nq. As before, I\ < For an upper bound on I 2 , we again apply Lemma 6.1 
in [3| to obtain 


h = 


OO 

£ 

n=N-\-l 


{fa fb’> u n) u n\\L 2 (T) 1 


<JN 


T 2 - 
1 2 ~ 


00 

£ 

n=N +1 


\(fa ~ fb,U n )\ 2 < C 2 (2k ) 2 jr 


This yields 


II fa — fb\\L 2 (F) < 2(5A 


n=N +1 


~ l e aN + 2 kc 


n z TV 


Vn' 


The value N = N(5) that minimizes the right-hand side in the above satisfies 

26A- 1 ae aN ^ - kcN(5 )~ 3/2 = 0 
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and therefore 

aN(S) + ^ log N(5) = log (^) • 
As a consequence, the optimal N ( 5 ) satisfies 


(a + |)JV(<y)>log(^) >aN(6) 


and hence, 


Thus, 


N(8y 


-I , 3, 

2 <{a + 


log 



1 

2 


(4.4) \\f a -f b \\ L 2 {F) < ^N(5)- 3 2+2kcN(5)- 1 i < kc(^ + 2 )N(5)~ 1 2 < K C[log(^) + D] * 

for some constants C > 0 and D that depend on 01 , 02 , 03 , 04 . This bound holds for 
N(5) > No, for which 

l < ^ e -(a+l)JV 0 

k 2 a 

is a sufficient condition. Note that by this requirement, log(^) + D > 0. 


5. Numerical validation of the asymptotic expressions in (2.5) 


We have calculated the SVD of the truncated Hilbert transform with oi = 0, 02 = 
450, 03 = 1350, 04 = 1725. The problem has been discretized with a sampling step equal 
to 1 and a shift of 1/2 between the object and data samples. Mathematica (using 20 
significant digits) finds 910 nonzero singular values for the discretized 1351 x 1276 Hilbert 
matrix Hx- As observed in [Ij, most singular values are very close to 1 , and only 10 
singular values are smaller than 0.97; of those the last nine are smaller than 0.01. These 
last singular values are well fitted by the asymptotic expression Ae~ an (Figure 2a) with 


the constants from equation (2.3). 

For each singular function u n corresponding to these last 9 values, we have computed 
the norm Wx^uWl 2 ^) i n th e interval ( 02,03 — //), for integer values of fi between 1 and 
400. For all fi, we find that |X/u' u n|| L' 2 (F) is we ll fitted by the asymptotic expression (2.5) 
as illustrated in Figure 2b. The fits are obtained by associating n = 1 to the singular value 


n 


= 902 of the discrete SVD. 
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Figure 2. Left: iogo+ versus n for the discretized truncated Hilbert problem (dots) 
with ai = 0, a2 = 450, a 3 = 1350 ,04 = 1725. The red line is the asymptotic expression 
(2.2 I. Right log ||XM M ™|ln 2 ([a 2 ,a 4 ]) versus n for g = 5 (red), /x = 20 (blue) and /x = 100 
(green) . The lines are the corresponding asymptotic expressions ( |2.5| ). The plots show 
the last nine singular components. 


6 . Appendix 


This appendix gives the proof of equations (2.5) and (2.6), which characterize the as¬ 
ymptotic behavior as n -> 00 of the integral 

/ ra 3—/x \ l/ 2 

(6-1) \\Xn u n\\L 2 (F) = ( / dx\u n (x) I 2 


where u n is the n-th singular function of the truncated Hilbert transform with overlap Ht- 
The proof relies heavily on the results from mm- These findings rely on the fact that 
the operator Ht shares an intertwining property with second-order differential operators 
Ls and L$ [2J. It implies that the singular functions u n and v n are specially constructed 
solutions to 


( 6 . 2 ) P(x)f"(x) + P'{x)f'{x) = (A n - 2(x - a) 2 ) f(x) 

where A„ -> 00 as n -> 00 and A n —> —00 as n —> —00 (see hkbi for details). Below we 
omit detailed references to these papers, to which the reader is referred for all properties 
stated without explanation. 


Lemma 1. For small £ > 0 the asymptotic behavior as n —> 00 of the singular function in 
x £ (02 + 0(e 1+2 ^),<23 — 0(e 1+21 ’)) is 


(6.3) 


U n (X = 


2 (- 1 ) 


71+1 


K. (P(x)) 1 / 4 


e -w 3 (x)/s( 1 + 0( e l/2-C)) 



dt 

VW) 


(6.4) 


with ws(x) 
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and where e is related to n by 
(6.5) 


K. 


e = - (1 + 0 (n _1/2_K )). 


T17T 


Proof. In the interval (ai + 0(e 1+2 ^), a 2 — 0(e 1+2 ^)), where e is related to n by (6.5), the 
WKB approximation of the singular function v n can be written as (see eqns (5.12) and 
( 8 . 1 ) in m 

Vn(x) = TST J£phw [ { le ‘ W,(X>, ’~‘ 1 ' ,a + e _ ‘" I< ' l/ ' + " /4 ](i + 0(^ 1/2 - c )) 

(6.6) + (!/j)[ e ”iW/ £ - CT / 4 - e -wiWA-HV4j 0 ( e i/2-C)J 


with 



dt 


Let v be the analytic continuation of v n from the interval (ai, 02 ) to C\[a 2 , 04 ]. Employing 
the uniqueness of the solution to the Riemann-Hilbert problem and the Plemelj-Sokhotski 
formulas, one can show (cf. eqns (4.2)-(4.8) in [1]) that 


, . v(x + *0) - v(x — iO) , , , 

u n (x) = a n -—-, x 6 (a 2 , a 4 )\{a 3 }, 

Zi 

. . v(x + zO) + v(x — *0) , , , 

V n (x) = ---, X € (ai,a 3 )\{a2j. 


Furthermore, we have that (cf. eqns (4.9), (4,10) in [1]) 


(6.7) u n (x) = cr n Im v(x + *0), r G (a 2 ,a 4 )\{a 3 }, 

( 6 . 8 ) v n (x) = Re v(x + *0), x € (ai,a 3 )\{o 2 }. 


The analytic continuation of —P(x) from ( 01 , 02 ) to ( 02 , a 3 ) via the upper half plane is 
given by —P(x) = e~ l7r P(x). With this, the analytic continuation of v via the upper half 
plane to (02 + 0(e 1+2< ’), a 3 — //) C (02 + 0(e 1+2< =), a 3 — 0(e 1+2< ’)) can be expressed as (see 
section 5.1 and Figure 2 in [lj, and m for the uniform accuracy of this continuation): 


v(x + iO) 


__ _ „»r/4 f \ p iK-/e—i-K/A ie x ' K / 2 'W 2 {x)/e 

V2 Kl(P(x)y/* V 

■ (1 + 0{eV 2 -t)) 


_|_ e —iK-/e+iir/A (D —ie l ' n t 2 W 2 {x)/£ 


(6.9) 


+ 


(1 /j)[ e *- ft: -/ e_ * 7r / 4 e *e l7r/2 «J2(a:)/£ _ e -iK-/£+iix/A e -ie™/ 2 W2{x)/£^Q^ £ A/2-C,^ 


with 


w 2 (x) = K + - w 3 (x) 



dt. 

7 wy 
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Equation (6.9) can be rewritten as: 

V ( T _i_ ;ni — - _ - f \p iK -l e p~ w ^ x )l e p -iK-/e+in/2 w 2 (x)/e l 

V[x + M) ~ V2K^(p(x)y/*i [ + 1 


( 6 . 10 ) 


• (1 + 0 (e 1/2 " c )) + (1 /i)[e iK - /£ e~ W2{x)/£ - e -i*-MW 2 e «*(*)/e] 0 ( e i/ 2 -C) J 


Noting that the modulus of all complex exponenti als is equal to 1 and that p w ‘F x )/ e > 
e —w 2 (x)/e^ one can regroup all terms in 0 {£ l ^ 2 ~^) in ( 6 . 10 ) to: 

(Q ll) _ - _ e W2(a;)/ir O('e 1/2_c ) 

[ ’ (P(x))V 4 1 J 

^ e ~iK~/e+i-K /2 e w 2 {x)/e _|_ & iK _/e e ~w 2 (x)/e + e w ^ x)/£ 0 {£ 1/2 -t)) 


Therefore, (6.10) becomes 

1 


v(x + iO) 

( 6 . 12 ) 


1 


V ^ K ^)) 1 / 4 

1 1 


y/2Kl (P(x )) 1 / 4 
Combining (6.7) and ( 6 . 12|) one obtains, 


e w 2 (x)/e ^ e -iK_/e+in/2 (^(g 1 / 2 -^). 


(6.13) 


(Jf, 


(_ i \n+l 


where we have used cos (K_/e) = (—l) n+1 (l + ©(e 1 / 2- ^)) (see eq. (5.36) in [I]). Finally, 
using the explic it ex pression (see (5.45 ) in (TJ) for the asymptotics of the singular values, 
yields equation (6.3) with w^(x) as in (6.4). This result holds for x E (02 + 0(£ 1+2(, ),O3 — 

0(e 1+2C )). □ 

Lemma 2. For sufficiently large positive index n, the singular function u n is strictly mono¬ 
tonic on ( 02 , 03 ). 


Proof. The singular function u n {x) is a solution of ( 6 . 2 ) for x E (02,04)^03}, where 
\ n — > 00 as n -> 00. By definition, P{x) > 0 for x E ( 02 , 03) and P'(a, 2 ) > 0. Since u n (x) 
is bounded at x = 02 and u n (a 2 ) f 0 (see e.g. | 2 j, Sec. 2 ), we may assume without loss of 
generality that u n {a 2 ) > 0. Then if 


(6.14) 


X n > max 2 (x — cr)* 

CL2<X<CL3 


the RHS of ( 6 . 2 ) is positive at 02 and since P{af) = 0, u'^af) is finite and P'{a 2 ) > 0, one 
concludes that u' n (a 2 ) > 0. We will now show that u' n {x ) > 0 for x E (02, 03). 

Let x be the first point after 02 for which u' n {x) = 0. Then, by construction we have 
u n {x) > 0 and u'f{x) < 0. With this, evaluating ( |6.2[ ) at x = x results in a contradiction 
since the LHS is non-positive whereas the RHS is positive. Hence, u' n (x) > 0 for all 
x E (02,03), i.e., u n is strictly monotonic. □ 


Proof of ( fg.5| ), pT6| ). Split the integral ( |6.1[ ) into 

r a 2+p 

dx\u n (x)\ 2 + / dx\u n (x)\‘ 
J aot 


(6.15) 


\\XpUr 


I L*{T) 


— la + h — 


raz-fL 


CL2+P 
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with p = 0 ( £ 1+2< >) and p, > 0 such that the integration interval for I a is contained in the 
domain of validity of the WKB approximation for n > N^. Using Lemma [lj the first term 
is 


Ia = 


f-a 3 —p 

' CL2+P 
£ 


dx \u n {x )\ 2 = ^( e - 2 ^ 3 (a 3 -M)/ £ _ e - 2 w 3 (a 2 +p )/ £ ) (1 + Q^ 1 / 2 ^)) 


K_ 


— c -2iD 3 (ffl 3 -/j)/e^ _ e 2/£[w 3 (a 3 -p,)-w 3 (a2+p)\^^ _|_ 


K_ 

£ 


= —e-^as-ri/eQ + 0 ( £ i/ 2 -C)), 
K 


where in the second line the first factor in parantheses gets absorbed in the error term since 
^ 3(^3 — pt) — 1 ^ 3(02 + p ) is negative and strictly bounded away from 0. Consider now the 
second term Iy. For sufficiently small e = 1 /\ /r X n , u n is monotonic on ( 02 , 03 ) by Lemma 
[2] Thus, 

ra 2 +p ra 3 -p, 

(6.16) h= dx\u n (x)\ 2 < - / dx\u n (x)\ 2 < 0(e 1+2c, )/ a 

J CL 2 ^3 fA ^2 p JCL 2 -\-p 

so that the contribution of can be absorbed in the error term of equation ( |3.5| ). We 
finally obtain 

/ - \l/2 

\\XpUn\\ L ^) = (la + h ) l/2 = (_ e - 2 -3(«3 -^)/ 6(1 + 0( £ V2-C )} j 

(6.17) = (1 + 0(n“ 1/2+c )) 

V n7r 

□ 
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